4  Regresion Lineal Simple

5 Introducción

Este documento presenta una guía detallada para el análisis de Regresión Lineal Simple a través de dos ejercicios prácticos. Se abordan todos los pasos del modelado: estimación de parámetros, pruebas de hipótesis con sus sistemas explícitos, construcción de intervalos de confianza y predicción, y la validación rigurosa de los supuestos del modelo.

5.1 Librerías Requeridas

library(MASS)
library(lmtest)
library(tseries)
library(ggplot2)

---

6 EJERCICIO 1: Relación entre Frecuencia Cardiaca y Peso

Se explora la relación entre la frecuencia cardiaca en reposo (frecuencia) y el peso corporal (peso).

6.0.1 Carga y Exploración de Datos

pacientes <- read.csv("frecuencia.csv", sep=";")
frecuencia <- pacientes$Frecuencia
peso <- pacientes$Peso
pacientes <- data.frame(frecuencia, peso) 
head(pacientes)

  frecuencia peso
1        127   75
2        101   46
3        127   79
4        143   97
5        104   48
6        107   46

plot(pacientes$peso, pacientes$frecuencia, xlab="Peso", ylab="Frecuencia", main="Diagrama de Dispersión")

6.1 A) Modelo Propuesto de Regresión

Se ajusta un modelo de la forma: \[ \text{frecuencia}_i = \beta_0 + \beta_1 \cdot \text{peso}_i + \epsilon_i \]

fit <- lm(frecuencia ~ 1 + peso, data = pacientes)
summary(fit)

Call:
lm(formula = frecuencia ~ 1 + peso, data = pacientes)

Residuals:
     Min       1Q   Median       3Q      Max 
-13.9734  -2.5950   0.8194   3.0306  13.3532 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 68.24820    1.98810   34.33   <2e-16 ***
peso         0.76295    0.02616   29.16   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 4.86 on 97 degrees of freedom
Multiple R-squared:  0.8976,    Adjusted R-squared:  0.8966 
F-statistic: 850.3 on 1 and 97 DF,  p-value: < 2.2e-16

Lectura del summary(fit): - Intercept (\(\hat{\beta}_0\)): 68.248 es la estimación del intercepto. - peso (\(\hat{\beta}_1\)): 0.763 es la estimación de la pendiente. - Std. Error para peso: Es la raíz de la varianza estimada de \(\hat{\beta}_1\). - Residual standard error: 4.86 es la estimación de la desviación estándar de los errores (\(\hat{\sigma}\)). - Adjusted R-squared: 0.8966 es el coeficiente de determinación ajustado. - degrees of freedom: 97 son los grados de libertad para los intervalos de confianza y pruebas de hipótesis de los coeficientes.

6.2 B) Estimación de Parámetros (Manual)

Se calculan manualmente los estimadores de los coeficientes, sus varianzas y la varianza del error.

\[ \hat{\beta}_1 = r_{xy} \frac{S_y}{S_x} \quad | \quad \hat{\beta}_0 = \bar{y} - \hat{\beta}_1 \bar{x} \quad | \quad \hat{\sigma}^2 = \frac{SCE}{n-2} \]

cor <- cor(pacientes$frecuencia, pacientes$peso)
s_y <- sd(pacientes$frecuencia)
s_x <- sd(pacientes$peso)
n <- length(pacientes$frecuencia)

b1_hat <- cor * s_y / s_x
b0_hat <- mean(pacientes$frecuencia) - b1_hat * mean(pacientes$peso)
y_hat <- b0_hat + b1_hat * pacientes$peso
sigma2_hat <- sum((pacientes$frecuencia - y_hat)^2) / (n - 2)
sigma_hat <- sqrt(sigma2_hat)

var_hat_b1_hat <- sigma2_hat / ((n - 1) * var(pacientes$peso))
var_hat_b0_hat <- sigma2_hat * sum(pacientes$peso^2) / (n * ((n - 1) * var(pacientes$peso)))
cov_hat_b0_b1 <- -sigma2_hat * mean(pacientes$peso) / ((n - 1) * var(pacientes$peso))
cor_hat_b0_b1 <- cov_hat_b0_b1 / sqrt(var_hat_b0_hat * var_hat_b1_hat)

# La función vcov() calcula la matriz de varianzas-covarianzas de los estimadores
vcov(fit)

            (Intercept)         peso
(Intercept)  3.95253557 -0.050422121
peso        -0.05042212  0.000684557

6.3 C) Intervalos de Confianza para los Parámetros \(\beta_0\) y \(\beta_1\)

Se construyen intervalos de confianza al 95% para los parámetros de localización del modelo. La fórmula es: \[ \hat{\beta}_j \pm t_{(1-\alpha/2, n-2)} \cdot EE(\hat{\beta}_j) \]

alpha <- 0.05
Li <- coef(fit) - qt(1 - alpha/2, n - 2) * sqrt(diag(vcov(fit)))
Ls <- coef(fit) + qt(1 - alpha/2, n - 2) * sqrt(diag(vcov(fit)))
cbind(Li, Ls)

                    Li        Ls
(Intercept) 64.3023726 72.194023
peso         0.7110212  0.814878

# Usando la función directa
confint(object = fit, level = 0.95)

                 2.5 %    97.5 %
(Intercept) 64.3023726 72.194023
peso         0.7110212  0.814878

6.4 D) Prueba de Hipótesis: ¿Hay relación entre las variables?

Se realiza una prueba de hipótesis para determinar si la pendiente \(\beta_1\) es significativamente diferente de cero.

Sistema de Hipótesis: \[ H_0: \beta_1 = 0 \quad \text{vs.} \quad H_1: \beta_1 \neq 0 \]

sign <- 0.05
t_c <- (as.numeric(fit$coeff - 0) / sqrt(vcov(fit)))

Warning in sqrt(vcov(fit)): NaNs produced

valcrit <- qt(df = n - 2, 1 - sign/2)

# Decisión por valor crítico
abs(t_c) > valcrit

            (Intercept) peso
(Intercept)        TRUE   NA
peso                 NA TRUE

# Decisión por p-valor
pval <- 2 * (1 - pt(abs(t_c), df = n - 2))
pval

            (Intercept) peso
(Intercept)           0  NaN
peso                NaN    0

Conclusión: Dado que abs(t_c) > valcrit es TRUE y el p-valor es prácticamente cero, se rechaza la hipótesis nula. Existe una relación lineal estadísticamente significativa entre el peso y la frecuencia cardiaca.

6.5 E) Prueba de Hipótesis: ¿Es la pendiente mayor que 1?

Sistema de Hipótesis: \[ H_0: \beta_1 = 1 \quad \text{vs.} \quad H_1: \beta_1 > 1 \]

sign <- 0.05
t_c <- as.numeric((fit$coeff - 1) / sqrt(vcov(fit)))

Warning in sqrt(vcov(fit)): NaNs produced

valcrit <- qt(df = n - 2, 1 - sign)

# Decisión por valor crítico
t_c > valcrit

[1]  TRUE    NA    NA FALSE

# Decisión por p-valor
pval <- 1 - pt(t_c, df = n - 2)
pval

[1]   0 NaN NaN   1

Conclusión: t_c > valcrit es FALSE y el p-valor es muy alto (0.999). No se rechaza la hipótesis nula. No hay evidencia para afirmar que la pendiente sea mayor que 1.

6.6 F) Intervalo de Confianza del 90% para la Media Condicional

Se estima el valor esperado de la frecuencia cardiaca para un peso de 70 kg.

alpha <- 0.1
x_new <- 70

mu <- coef(fit) + coef(fit) * x_new
var <- vcov(fit) + x_new^2 * vcov(fit) + 2 * x_new * vcov(fit)

Li <- mu - qt(1 - alpha/2, n - 2) * sqrt(var)

Warning in sqrt(var): NaNs produced

Ls <- mu + qt(1 - alpha/2, n - 2) * sqrt(var)

Warning in sqrt(var): NaNs produced

cbind(Li, Ls)

            (Intercept)    peso (Intercept)     peso
(Intercept)    4611.204     NaN     5080.04      NaN
peso                NaN 51.0844         NaN 57.25444

El intervalo resultante indica, con un 90% de confianza, entre qué valores se encuentra la frecuencia cardiaca promedio para todas las personas que pesan 70 kg.

6.7 G) Intervalo de Predicción del 99%

Un intervalo de predicción estima el rango en el que se espera que caiga una nueva observación individual.

alpha <- 0.01
x_new <- 85

sigma2 <- sum(fit$residuals^2) / (n - 2)
mu_pred <- coef(fit) + coef(fit) * x_new
var_pred <- sigma2 * (1 + 1/n + (x_new - mean(peso))^2 / ((n-1)*var(peso)))
Li_pred <- mu_pred - qt(1 - alpha/2, n - 2) * sqrt(var_pred)
Ls_pred <- mu_pred + qt(1 - alpha/2, n - 2) * sqrt(var_pred)
cbind(Li_pred, Ls_pred)

               Li_pred    Ls_pred
(Intercept) 5856.48662 5882.20340
peso          52.75527   78.47206

6.8 H) Descomposición de la Varianza (Tabla ANOVA)

La variabilidad total (SCT) se descompone en la explicada por el modelo (SCM) y la no explicada (SCE).

\[ \sum(y_i - \bar{y})^2 = \sum(\hat{y}_i - \bar{y})^2 + \sum(y_i - \hat{y}_i)^2 \]

SC_Tot <- sum((pacientes$frecuencia - mean(pacientes$frecuencia))^2)
SC_residuales <- sum(fit$residuals^2)
SC_mod <- sum((y_hat - mean(pacientes$frecuencia))^2)

all.equal(SC_Tot, SC_mod + SC_residuales)

[1] TRUE

6.9 I) Significancia de la Regresión y Coeficiente de Determinación

Se realiza una prueba F para la significancia global del modelo.

CM_residuales <- SC_residuales / (n - 2)
CM_mod <- SC_mod / 1
F_c <- CM_mod / CM_residuales

# Decisión por p-valor
p_valor_f <- 1 - pf(F_c, df1 = 1, df2 = n - 2)
p_valor_f

[1] 0

6.10 J) Validación del Modelo

Se verifican los supuestos del modelo sobre los residuales.

residuals <- fit$residuals

# a. Media cero y no relación con variables
par(mfrow=c(1,2))
plot(frecuencia, residuals, main="Residuales vs. Y")
plot(peso, residuals, main="Residuales vs. X")

mean(residuals)

[1] -3.385339e-17

t.test(residuals, mu = 0)

    One Sample t-test

data:  residuals
t = -6.9659e-17, df = 98, p-value = 1
alternative hypothesis: true mean is not equal to 0
95 percent confidence interval:
 -0.9644195  0.9644195
sample estimates:
    mean of x 
-3.385339e-17 

# b. Varianza constante (Homocedasticidad)
plot(y_hat, residuals, main="Residuales vs. Ajustados")
bptest(frecuencia ~ 1 + peso, data=pacientes)

    studentized Breusch-Pagan test

data:  frecuencia ~ 1 + peso
BP = 0.29529, df = 1, p-value = 0.5868

# El p-valor es grande, así que no se rechaza que haya homocedasticidad.

# d. Normalidad de los residuales
par(mfrow=c(1,2))

hist(residuals, main="Histograma de Residuales")
qqnorm(residuals); qqline(residuals, col="red")

jarque.bera.test(residuals)

    Jarque Bera Test

data:  residuals
X-squared = 2.1418, df = 2, p-value = 0.3427

shapiro.test(residuals)

    Shapiro-Wilk normality test

data:  residuals
W = 0.97912, p-value = 0.1173

# Los p-valores son grandes, por lo cual no se rechaza la normalidad.

Conclusión: Todos los supuestos del modelo parecen cumplirse satisfactoriamente.

---

7 EJERCICIO 2: Efecto de un Tratamiento para la Hipertensión

7.0.1 Carga y Preparación de Datos

hipertension <- read.table("hipertension.csv", header=T, sep=";")
diferencia <- hipertension$Antes - hipertension$Despues
tratamiento <- hipertension$Tratamiento
hipertension <- data.frame(tratamiento, diferencia)
head(hipertension)

  tratamiento diferencia
1        1PLA        1.0
2        1PLA        0.6
3        1PLA       -0.9
4        1PLA       -0.5
5        1PLA       -3.3
6        2MED       11.0

7.1 B) Estimación de Parámetros del Modelo

El modelo es: \(\text{diferencia}_i = \beta_0 + \beta_1 \cdot \text{tratamiento}_i + \epsilon_i\), donde tratamiento es 1 si es “2MED” y 0 si es “PLACEBO”.

fit2 <- lm(diferencia ~ 1 + tratamiento, data = hipertension)
summary(fit2)

Call:
lm(formula = diferencia ~ 1 + tratamiento, data = hipertension)

Residuals:
   Min     1Q Median     3Q    Max 
-3.420 -0.310  0.300  1.265  1.980 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)      -0.6200     0.8511  -0.728    0.487    
tratamiento2MED  11.1400     1.2037   9.255 1.51e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.903 on 8 degrees of freedom
Multiple R-squared:  0.9146,    Adjusted R-squared:  0.9039 
F-statistic: 85.66 on 1 and 8 DF,  p-value: 1.508e-05

7.2 C) Intervalos de Confianza del 99% para los Parámetros

alpha <- 0.01
n2 <- nrow(hipertension)
Li <- coef(fit2) - qt(1 - alpha/2, n2 - 2) * sqrt(diag(vcov(fit2)))
Ls <- coef(fit2) + qt(1 - alpha/2, n2 - 2) * sqrt(diag(vcov(fit2)))
cbind(Li, Ls)

                       Li        Ls
(Intercept)     -3.475827  2.235827
tratamiento2MED  7.101251 15.178749

7.3 E) Prueba de Hipótesis: \(H_0: \beta_1 = 0\) vs. \(H_1: \beta_1 > 0\)

Se prueba si el medicamento tiene un efecto significativamente mayor que el placebo.

t_c <- (as.numeric(fit2$coeff - 0) / sqrt(vcov(fit2)))

Warning in sqrt(vcov(fit2)): NaNs produced

pval <- 1 - pt(t_c, df = n2 - 2)
pval

                (Intercept) tratamiento2MED
(Intercept)       0.7564458             NaN
tratamiento2MED         NaN    7.541131e-06

7.4 G) Comparación con un Modelo ANOVA

Una regresión con una sola variable predictora categórica es matemáticamente equivalente a un ANOVA de una vía.

fit_aov <- aov(diferencia ~ tratamiento, data = hipertension)
summary(fit_aov)

            Df Sum Sq Mean Sq F value   Pr(>F)    
tratamiento  1 310.25  310.25   85.66 1.51e-05 ***
Residuals    8  28.98    3.62                     
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

summary(fit2)

Call:
lm(formula = diferencia ~ 1 + tratamiento, data = hipertension)

Residuals:
   Min     1Q Median     3Q    Max 
-3.420 -0.310  0.300  1.265  1.980 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)      -0.6200     0.8511  -0.728    0.487    
tratamiento2MED  11.1400     1.2037   9.255 1.51e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.903 on 8 degrees of freedom
Multiple R-squared:  0.9146,    Adjusted R-squared:  0.9039 
F-statistic: 85.66 on 1 and 8 DF,  p-value: 1.508e-05